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ABSTRACT 


An abstract approximation framework is developed for the finite and 
infinite time horizon discrete-time linear-quadratic regulator problem for 
systems whose state dynamics are described by a linear semigroup of operators 
on an infinite dimensional Hilbert space. The schemes included in the 
framework yield finite dimensional approximations to the linear state feedback 
gains which determine the optimal control law. Convergence arguments are 
given. Examples involving hereditary and parabolic systems and the vibration 
of a flexible beam are considered. Spline-based finite element schemes for 
these classes of problems, together with numerical results, are presented and 
discussed. 
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1. Introduction 


Recent advances in micro-processor technology have led to increased 
interest in digital or discrete-time control systems. In addition, because 
many current application areas involve complex systems which are most 
appropriately modelled using functional and/or partial differential equations, 
it has become important to study digital control techniques in the context of 
infinite dimensional or distributed systems. 

A great deal of attention has been given to the continuous-time infinite- 
dimensional linear-quadratic regulator problem. The general theory and 
characterization of the linear state feedback form of the optimal control are 
discussed in [5], [6], [8], [9], [21] and [22], while its application to 
hereditary, parabolic and hyperbolic systems with emphasis on approximation is 
treated in [2], [3], [7], [10], [11], [14] and [17] to mention just some of 
the work that has been done. 

On the other hand, relatively little can be found in the literature 
concerning the corresponding discrete-time problem. The major contributions 
in this area can be found in the papers by Lee, Chow and Barr [20] and Zabczyk 
[28]. In these studies the Riccati difference equations that characterize the 
linear feedback form of the optimal control for the finite time problem are 
given and limiting properties as the length of the time horizon tends to infinity 
are discussed. However, the issue of approximation is not considered. 

In the present paper, we develop numerical approximation schemes that 
yield finite dimensional approximations to the feedback gain operators which 
determine the discrete— time optimal control law. We consider control systems 
whose dynamics can be described in terms of a linear semigroup of operators on 
an infinite dimensional Hilbert space. The basis of our approach is the 
construction of a sequence of finite dimensional (presumably finite element 
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based) state approximations which in turn leads to a sequence of finite 
dimensional discrete-time linear-quadratic regulator problems each of which 
can be solved using standard techniques. 

Under appropriate assumptions on the nature of the original problem and 
the convergence of the state approximation, we are able to prove that the 
approximating optimal controls and feedback gains converge to the true optimal 
control sequences and feedback laws for the original infinite dimensional 
system. Depending upon the convergence properties of the state approximation, 
we are able to establish strong or uniform norm convergence of the 
approximating gain operators and the corresponding weak or strong convergence 
of the approximating feedback kernals which are used in the implementation of 
the optimal control. We treat both the finite and infinite-time horizon 
problems. 

We have tested our schemes on a wide variety of examples. This paper 
includes numerical results for problems with state dynamics given by 
hereditary and parabolic (heat/diffusion) differential equations and a hybrid 
system of partial and ordinary differential equations for the vibration of an 
Euler-Bernoulli beam connected to a rigid body and a lumped mass. We 
implemented and tested the methods on an IBM Personal Computer. 

We give a brief outline of the remainder of the paper. In section 2 we 
breifly outline previous results concerning the characterization of the 
optimal control and feedback gains for both the finite and infinite time 
horizon discrete-time regulator problem for distributed systems. The Riccati 
difference and algebraic equations whose solutions determine the optimal 
feedback control law are discussed. In section 3 we develop the abstract 
approximation framework and convergence arguments. Section 4 contains a 
discussion of particular schemes for the classes of problems mentioned above 
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together with the results of our numerical studies. Some concluding remarks 
are given in Section 5. 

We employ standard notation throughout. For an interval (a,b), we denote 
by H (a,b) the usual Sobolev spaces of real-valued functions defined on (a,b) 
whose (k-l)st derivatives are absolutely continuous and whose k 1 "* 1 derivatives 
are L^. The standard Sobolev inner product on H (a,b) is denoted 
by <*,*>^. For X and Y normed linear spaces we denote by 1(X,Y) the space of 
bounded linear operators from X into Y. When Y = X, we use the shorthand 
notation L(X). 

2. The Optimal Control Problem 

2. 1 Optimal Control on a Finite Interval 

Let Z and U be Hilbert spaces with inner products <•,•> and <•,•>.., 

u 

respectively, with U finite dimensional. For {H, <•,•>} a Hilbert space, let 

2 t £ 

1 (tQ,t^;H) denote the usual Hilbert space of sequences x = {x(t)} with 

r t t Q 

x(t) e H together with the inner product 

fc f 

(2.1) <x,y> _ = l <x( t ) , y(t)> . 

* t-c o 

The discrete-time linear quadratic regulator problem on the finite time 
interval [tQ,t^] is 

2 

(Pi) Choose u g % (tQ,t^;U) to minimize the quadratic performance index 
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J ^ G;t 0 ,t: f » z ^ t o^ ,u ^ = 

(2.2) t f -l 

I [<Qz(t) , z(t)> z + <Ru(t) ,u(t)>y ] + <Gz(t f ),z(t f )> z 
t=t 0 

subject to the discrete-time control system 


z(t+l) = Tz(t) + Bu(t), t > t 

(2.3) 

z(t Q ) e Z, 


where T and B are bounded linear operators from Z into Z and U into Z, 
respectively, Q and G are bounded, nonnegative self-adjoint operators on Z, 
and R is a positive definite self-adjoint operator on U. 

Of primary concern to us will be applications where (2.3) is the sampled 
form of the continuous-time control system 

(2.4) z(s) =Az(s) + Bu(s) 

where A is the infinitesimal generator of a Cg-semigroup of bounded linear 
operators T( s), s > 0, on Z, and B is a possibly unbounded linear operator 
from U into Z. In this case we have 

(2.5) T = T(x ) and B = /J T(s)Bds, 

where t is the sampling interval. If, as in our subsequent example 
discussed in Section 4.1 where u is a boundary control in a heat equation, B 
is unbounded (more precisely, B maps U not into Z but into some larger space), 
then the integral in (2.5) is not interpreted literally. 
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The solution to Problem (Pi) has been given for infinite dimensional 

control systems in [20], [28], and the equations representing the solution have 

the same form as in the finite dimensional case. We will give now the version 

of the solution that is most useful for our purposes. 

For given z(tg), J(G;tg,t^,z(tg),u) is a bounded linear-quadratic 
2 

functional on Z (tg,t^;U) with coercive quadratic part. Therefore, for each 

2 

z ( tg ) , there exists a unique optimal control sequence in Z (tg,t^;U). Also, 
the minimum value of the performance index is a quadratic functional of z(tg), 
so that there exists a unique nonnegative, self-adjoint n(tg) e L( Z) such 
that 

(2.6) J* = min J(G;t Q ,t f ,z(t Q ),u) = <n(tg)z(tg),z(tg)> z . 

Application of the principle of dynamic optimality establishes that the 
optimal control has the feedback form 

(2.7) u*(t) = -F(t)z*(t), tg < t < t f -l 

where 

(2.8) F(t) = R(t) _1 B*n(t+l)T, 

(2.9) R(t) = R + B*n(t+1)B 

and n(t) satisfies the Riccati difference equation 

(2.10) n(t) = T*[n(t+i) - n(t+i)BR(t) - 1 B*n(t+i)]T + q. 


t < t f -l. 
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with the final condition 


(2.11) n(t f ) = G. 


The optimal trajectory z* is given by 


(2.12) z*(t+l) = S(t)z*(t), t > t Q 


where 

(2.13) S(t) = T-BF(t). 

2.2 Control on the Infinite Interval 

Here, t^ = 00 and G = 0. To simplify notation, we will write 
J ( t O» 00 » z ( t o)» u ) instead of J(0; t Q ,“,z(t Q ),u). 

2 

Definition 2.1. A control sequence u e l (0,°°;U) is an admissible control 
for the initial condition z if J(0,°°,z,u) < 00 . 

The discrete-time linear-quadratic regulator problem on the infinite 
interval is 

(P2) Choose an admissible control u* to minimize J(0,°°,z,u), if an 

admissible control exists for the initial condition z. 

That a unique optimal control u* exists whenever at least one admissible 

control exists follows from the fact that the quadratic part of J(0,“,z,u) is 

2 

coercive on a subspace of l (0,°°;U). See the discussion following Definition 
4.1 of [9]. 
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Definitlon 2.2. A bounded linear operator II on Z is a solution to the 
Riccati algebraic equation if 

(2.14) n = T*[n-nB(R+B*IIB) -1 B*n]T + Q. 

The following theorem summarizes results from Zabczyk [28]. 

Theorem 2.3. The following are equivalent: 

(i) There exists an admissible control for each z e Z; 

(ii) for each z e Z, sup <II(t)z,z> < 00 , where II ( t ) is the Riccati 

C<t f 

operator in (2.10) and n(t^) = 0 for fixed t^; 

(iii) as t - 00 , n(t) converges strongly to a nonnegative self-adjoint 
solution to the Riccati algebraic equation; 

(iv) there exists a nonnegative self-adjoint solution to the Riccati 
algebraic equation. 

For uniqueness of the solution to the Riccati algebraic equation and 
characterization of the optimal control, Zabczyk treated two cases: when Q is 

coercive, and when the spectral radius of T is less than 1 (i.e. , the open- 
loop system is uniformly exponentially stable). Since neither is the case in 
the example we discuss in Section 4.2 and other applications in which we are 
interested, we will need the following hypothesis and theorem. 

Hypothesis 2.4. The operators T, B and Q are such that, if z(0) e Z and u is 
an admissible control for z(0), then 

lim |z(t) L = 0. 

f*CO 


(2.15) 
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Theorem 2.5. When Hypothesis 2.4 holds, there exists at most one nonnegative 
self-adjoint solution to the Riccati algebraic equation. If such a solution 
II exists, then there exists a unique solution to problem (P 2) for each 
initial condition z(0) e Z, the minimum value of the performance index is 

<2 - 16> J * ‘ » admissible J(0 -' ’ <H*«»,z<0» z , 
the optimal control has the feedback form 

(2.17) u*(t) = t > 0, 


where 

(2. 18) F = R _1 B*n T, 

(2.19) R = R+B*II B 

and the optimal trajectory z*(t) satisfies 

(2.20) z*(t+l) = Sz A (t), t > 0, 


with 


(2.21) S = T-BF. 


Proof. Let II be such a solution and note that, for any finite tf, II is a 
constant solution to (2.10) and (2.11) with G = II. Then the corresponding 
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F(t) and R(t) defined by (2.8) and (2.9) are the constant operators in (2.18) 
and (2.19). For z(0) £ Z, define z(0) = z(0), 

(2.22) 7(t+l) = (T-BF )z"( t ) , t > 0, 

and 


(2.23) u(t) = -Fz(t), t > 0. 


Now suppose that u is an admissible control for z(0) and that z(t) is the 
corresponding solution to (2.3). For t £ > 0, the preceding results about the 
solution to Problem (Pi) with G = II imply 

J(H;0,t f ,z(0),u) < J(0;0,t f ,z(0),u) + <IIz(t f ),z(t f )> z 

(2. 24) 

< J(0;0,°°,z(0) ,u) + <Hz(t f ),z(t f )> z . 


Also, 


(2.25) 


J(II;0,t f ,z(0),u) = <IIz(0),z(0)> z 


= J(0;0,t f ,z(0),u) + 


<01z(t f ),z(t f )> z . 


Since z(t^) + 0 as t^ + 00 , (2.24) shows that u is both admissible and 

optimal for Problem (P 2). Since z(t^) -*■ 0 as t^ + (2.25) shows (2.16). 

As we see now, (2.16) must hold for any nonnegative self-adjoint solution of 
the Riccati algebraic equation; therefore, such a solution is unique. 


Remark 2.6. When Hypothesis 2.4 does not hold, the Riccati algebraic equation 
may have more than one nonnegative self-adjoint solution. In this case, the 

minimal such solution — there will be one — gives the solution to Problem 
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(P2) as in Theorem 2.5. Throughout this paper, we assume that Hypothesis 2.4 
holds. 


Lemma 2.7. Suppose that Q > m for some positive constant m, and set 
n * 

C = Z (T ) QT , for n = 1,2, ... . Then |C z| 7 is bounded in n for each 
t=0 

z e Z if and only if C n converges in norm to the operator 

00 

(2.26) C = I (T*) t QT t 

t=0 


and 

(2.27) |T t | < ( | C | /m) (1 - m/lcl)*, t = 1,2 

Proof* Since C n is an increasing sequence of bounded self-adjoint linear 
operators, C n converges strongly to some bounded self-adjoint C if and only if 
<C n z,z> z is bounded in n for each z, if and only if |C n z| z is bounded in n for 
each z. This is a standard result. The proof of the Lemma is then a standard 
exercise using the Lyapunov functional <Cz(t), z(t)> z for the homogeneous part 
of (2.3). 

Corollary 2. 8. If Q > m > 0 and the Riccati algebraic equation has a 
nonnegative self-adjoint solution II, then the spectral radius of the operator 
S in (2.21) is less than 1, and 

(2.28) [S 11 ! < (|n|/m)(l - mini) 11 , t-1,2 


Proof. This follows from Lemma 2. 7 and 
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(2.29) <OIz,z>_ = Z (S*) t [Q + F^FlS 11 . 

L t-0 

For Q coercive, Zabczyk proved a stronger result than part (iii) of 
Theorem 2.3: if a nonnegative self-adjoint solution to the Riccati algebraic 

equation exists, then | II ( t ) - II | 0 geometrically fast as t + ■* (Also, see 

[13]). We will need such a result, along with an explicit convergence rate, 
for the approximation theory in Section 3.2. Since Zabczyk' s proof does not 
yield an explicit convergence rate, we give the following. 

Theorem 2.9. Suppose that there exists a nonnegative self-adjoint solution 
II to (2.14) and that 

(2.30) \s t \ < Mr t , t = 1,2, ... , 

where M and r are positive constants with r < 1 and S is the optimal closed- 
loop operator in Theorem 2.5. If II(») is the operator in (2.10) with t^ = 0 
and 

(2.31) n(0) > n. 


then 

(2.32) <IIz,z >2 < <n(-t)z,z> 2 < <IIz,z> z + (Mr t )^|n(0) | , t = 1,2, ... . 


Proof. For tQ a negative integer, let uq be the optimal control sequence for 
the finite-time Problem (Pi) on the interval [tg,0] with initial condition 
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z(t Q ) e Z, with Zq the corresponding optimal trajectory. Also, let u* be the 
optimal control sequence on the infinite interval for Problem (P 2) with 
initial condition z(tQ), with z A the corresponding optimal trajectory. 

Since II is a constant solution to (2.10) for the final condition G =11 , we 
have 


(2.33) 


< nz (t o ) , z ( t 0 ) > z = J(n ; t 0 ,o,z(t 0 ),u A ) 

< J(0;t 0 ,0,z(t 0 ),u 0 ) + <nz 0 (0),z 0 (0)> z 

< J(0;t 0 ,0,z(t 0 ),u 0 ) + <n(0)z 0 (0),z 0 (0)> z 


= <n(t 0 )z(t 0 ),z(t 0 )> 2 . 


On the other hand (note that z A (tg) = S z(tg)), 

<n(t 0 ) z (t 0 ),z(t 0 )> z 

< J(0;0,-t 0 ,z(t 0 ),u A ) + <n(0)z*(-t 0 ),z*(-t 0 )> z 

< J(0;0 > » > z(t 0 ),u*) + <n(0)z*(t 0 ),z*(t 0 )> z 

< <nz(t 0 ),z(t 0 )> z + | n (o) | ( | s to | |z(t 0 )| z ) 2 . 


(2.34) 
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3. Approximation Theory 


3.1 The finite time interval problem 

In this section we develop a general approximation framework for the 

finite time interval problem (Pi) and describe associated convergence results. 

For each N = 1,2, ... , let Z N c z be a finite dimensional subspace of Z 

and let P„: Z Z„ denote the orthogonal projection of Z onto Z M with respect 
N N 

to the <*,*>„ inner product. We require the following hypotheses. 

Cj 


Hypothesis 3.1 There exist operators T N : Z N Z N> B N : U + Z N 

Q. t : Z, T Z„ and G M : Z M -*• Z XT which satisfy 

N N N N N 


l N‘N 

* 

[ N P N 

3 n" 

^N P N 

3 n p n 


-»• T 

strongly. 

* 


r" T 

strongly, 

B 

strongly. 

♦ Q 

strongly. 

G 

strongly. 


as N 00 with T N and B N bounded and Q N and G N bounded, self-adjoint and 
nonnegative. 


Hypothesis 3.2 The spaces Z N are approximating subspaces in the sense that 
the projections P N satisfy P N -► 1 strongly on Z as N + ®. 

We note that since U has been assumed to be finite dimensional, 

* * 

Hypothesis 3. 1 above necessarily implies that B^ -*■ B and B^P^ B in the 
uniform norm topology on L(U,Z) and L(Z,U) respectively. 
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We define a sequence of approximating discrete-time linear quadratic 
regulator problems on the finite time interval [tg,t f ] as follows: 


< Pl H> 


N 2 

Find u A e SL (tQ,t^-l;U) which minimizes 


t f“ 1 

(3 - 1) J N (G N ;t 0 ,t f ,z(t 0 ),u) = l [<Q N z N (t),z N (t)> z + 

t-t 0 


<Ru(t),u(t)> u 3 + <V N (t f ), z N (t f )> 2 


subject to 

(3.2) z N (t+l) - T K z N (t) + BjjU(t) , t > t Q 


Z N <C 0^ ‘ V (t 0>- 


The results stated in Section 2. 1 concerning the existence and uniqueness of 

solutions to Problem (PI) apply to the Problems (Pl N ) as well. Indeed, there 

N 2 

exists a unique solution u* e Z (t Q ,t f -l;U) to Problem (P1 N ) which is given 
in feedback form by 

(3.3) uj(t) = -F N (t) zj(t), t Q < t < t f -l 


where 

(3.4) Fjj(t) “ R N (t) -1 B* Hjj(t+1)T H 
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with 

(3.5) R N (t) = R + B* n N (t+l)B N 

and the operators {II (t)} on Z N satisfying the Riccati difference 

N t t Q 

equation 

(3.6) n N (t) - yn„Ct+i) - n N (t+i)B N R N (t)' 1 B^n N (t+i)]T 1J + Q„ 

with terminal condition 


(3.7) n H (t f ) - g n . 


N 

The optimal trajectory z. is given by 


(3.8) 


z*(t+l) = S N (t)z*(t), t > t Q , 


N 


Z * ^0^ P N Z ^ t 0^ 


where 


(3.9) S N (t) = T N - B N F N (t), t > t Q . 


The operators {II (t)} are bounded, self-adjoint and nonnegative. 

N t t Q 

minimum value of the performance index (3.1) is given by 


(3.10) J* J^(G N ;tQ,t^,z(tQ),u ;fc ) <IlN(to)z*(to)» z *^ t 0^' > Z * 


The 


The fundamental convergence result is given in the following theorem, 
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N 

Theorem 3.3 Let u* and u* be the unique solutions to problems (Pl N ) and (PI), 

N 

respectively, with z* and z* the corresponding optimal trajectories 
generated by (3.8) and (2.12). Let Jjj, 11^ and F^ and J, II and F be given by 
(3.1), (3.6) and (3.4) and (2.2), (2.10) and (2.8). Then, if Hypotheses 3.1 
and 3.2 hold, we have 


(i) 

(ii) 

(iii) 

(iv) 


lim 

1 N 1 
l u * " U J 2 

= o. 

N-h» 

jr 


lim 

1 N 1 
l z * ~ z *l 2 

= o. 

N+oo 

JT 


lim 

14 - J*l = 

o. 




lim 

N+°° 

|n N (t) V - 

n(t)z| z 


and 


(v) lim | F (t )P - F(t) | =0, t n < t < t_ -1. 

w w u t 

Proof 


We first note that n^(t) being nonnegative implies that |R^(t)| > ( R | 
and consequently that I ^ N ( t) | < |r| . It follows therefore that for 

u e U. 


(3.11) |(R N (t) _1 - R^rSuly 


= 1^(0 1 (R(t ) - R N (t))R(t)" 1 u| u 
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< |R| _1 |(R(t) - R N (t))R(t) 

The above estimate together with (2.9), (2.11), (3.5), (3.7) and Hypothesis 
3.1 imply that 

(3.12) R N (t f -l) _1 + R(t f -1 ) _1 

as N ->■ " strongly on U. Since U is finite dimensional the convergence in 

(3.13) is in fact uniform. It then follows immediately from (2.8), (3.4) and 
Hypothesis 3.1 that 

(3.13) F N (t f - 1) P N - F(t f - 1), 
uniformly as N + ", and from (2.10) and (3.6) that 

(3.14) n N (t f - i)p n ♦ n(t f -i) 

strongly on Z as N + A simple induction yields (iv) and (v) from which 

(i), (ii) and (iii) then follow trivially. 

Remark It will, on occassion, be the case that in constructing a particular 

* * 

approximation scheme -*■ T strongly but T N P N *■ T only weakly (see, for 

example, [3]). However, by using the fact that 

( W t+1)) * "V t+1)T H 


(3.15) 
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^ & 

implies that T^Ct+l) Mfl (t+1) weakly if II^t+1) ^ n(t+l) weakly, we 
conclude that Theorem 3.3 continues to hold under these somewhat weaker 
hypotheses with the strong convergence in (iv) replaced by weak and the 
uniform convergence in (v) replaced by strong. 

Under certain additional hypotheses it can be shown that the operators 
n (t)» < t < tf given by (2.10), (2.11) are trace class (see [15]) and that 

(3.16) lim nn N (t)P N - n(t)n 1 = o, t Q < t < t f , 

where II • II ^ denotes the trace norm, the strongest of all common operator 

norms. We require the following lemmas. 

00 

Lemma 3.4 If is an absolutely summable sequence of real numbers then 

00 00 oo 

there exist sequences i i==1 and {c^^ such that lim b ± = 0, is 

absolutely summable and a^ = b^c^. 

Proof 

Let 

00 

(3.17) a = l | a | 

i=l 1 

and for j = 0,1,2,... define nonnegative integers nj as follows. Let nQ = 0 
and let nj denote the first index for which 

n . 

(3.18) |a | > a - \ , 

i=l i 

j ~ 1 • • • • Set 
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(3.19) b i = J » c i = J a i» 1 = n j-i +1 »*** » n j » j = 1 » 2 >' 


Then b.c. = a., i = 1,2,..., lim b. = 0 and 
l i i ’ . i 


(3.20) l |c | = l j l |a | < a + £ \ < - . 

i=l j=l k«n +1 K j=l 


Lemma 3.5 If L is a self-adjoint trace class operator on a separable Hilbert 

12 1 2 

space H, then L can be written as L L where L is compact and L is trace 
class. 


Proof 

°0 

Let denote the eigenvalues of L repeated according to 

00 

multiplicity and let {<|>^}^_^ denote the corresponding eigenvectors. Then 

CO 

is a sequence of real numbers, each of finite multiplicity, and 


(3.21) 


l l\l - Wj < 

i-i 


00 00 

Applying the previous lemma there exist sequences {y } , and {v } 

oo i 1 i. 1 1 X 

lim y . = 0, £ |v | < °° and X. = y.v . Defining L* and L^ by 

i+- i=l ill 

00 

(3.22) L 1 ^ = l » 1 i < *»*iVi» * e H 


(3.23) L 4> = l v <<(.,<(, > (j> <(> e H 

i=l 11 


with 
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respectively, the lemma immediately follows. 

00 

Lemma 3.6 Let t N=1 3 se< l uence °f bounded linear operators on a 

seperable Hilbert space H which converges strongly to a bounded linear 

oo 

operator S. Let be a sequence of trace class operators on H which 

converges in trace norm to an operator L. If L can be written as L = L^L^ 
with L compact and L trace class then the sequence { S N L N ^ N= ^ converges in 
trace norm to SL. 

Proof 

The result follows immediately from 

(3.24) II s n l n - SLIj < ns N (L N - Dlj + II (S N - SJL^Ij 

< l S N l' IL N " Lll l + l (S N " S )L 1 1 II L 2 II x . 

Theorem 3.7 If Q and G are trace class operators then the operators 

{n ( t )} given by (2.10) and (2.11) are trace class. Moreover, if 

0 

Hypotheses 3.1 and 3.2 hold and Q. T P. T Q and G P G in trace norm as 

N N N N 

N -*■ 00 then we have 

(3.25) lim nn (t)p - n ( t ) II =0, t n < t < t_. 


Proof 

That the operators II (t), tQ < t < t^ are trace class is an immediate 
consequence of the hypotheses of the theorem, (2.10), (2.11) and the fact that 
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the trace class operators form a two sided ideal of L(Z), the space of bounded 
linear operators on Z (see [15]). 

The trace norm convergence stated in (3.25) will follow once we have 
shown that 


(3.26) lim lin (t+l)P - n(t+l)ll 1 = 0 

N-h» 


implies 


(i) lim IITjl N (t+l)T N P N - T n(t+l)TII l = 0 
N+°° 

and 

(ii) lim IT* H (t+1) B N R N (t)" 1 B^I N (t+l)T N P N - 

T*n(t+l)BR(t) -1 B*n(t+l)TII = 0. 


To argue (i) we first note that Hypothesis 3.1 and Lemmas 3.5 and 3.6 


imply 


(3.27) lim IT* n„(t+l)P„ - T*n(t+l)l. = 0. 

N -H» N N N 1 


Taking ad joints we obtain 


(3.28) lim in„(t+l)T„P„ - n(t+l)TI. = 0. 

N-hx, N N N 1 

Another application of the previous two lemmas yields 


lim IT* II T (t+1)T„P„ - T*n(t+1)TI. < lim 
N N N N 1 


It I 

i N i 


in (t+i)T p - iKt+DTi. 
N N N 1 


(3.29) 
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+ lira |(T*P - T*) n 1 (t+i)|nn 2 (t+i)Tii. 

where n(t+l) = n 1 (t+l)n 2 (t+l) is the factorization of II(t+l) described in 
Lemma 3.6. 

The verification of (ii) is analogous and the theorem is proven. 

We note that if Hypotheses 3.1 and 3.2 hold and if the operators Q and G 
are trace class with Qjj and defined by 

(3.30) Q n = P n Q 
and 

(3.31) G n = P n G, 

then Lemmas 3.5 and 3.6 imply that the trace norm convergence hypotheses in 
Theorem 3.7 hold. The significance of this observation will become apparent 
when examples are discussed in Section 4. Indeed, it is frequently the case 
in practice that Q and G have finite rank (and consequently are trace class) 
and the operators Q N and G^ are defined as in (3.30), and (3.31). 

3.2 Approximation on the Infinite Interval 

Problem (P 2 N ) is Problem (P 2) for the control system in (3.2) and the 
performance index 


00 



-23- 


Hypothesis 3.8. For each N, there exists exactly one nonnegative self-adjoint 
solution to the Riccati algebraic equation 

(3.33) it N - yn N - Vn (r ♦ b b n )" 1 b> n ]t n + Q r 
By Theorem 2.3, this implies that 

(3.34) lim |n - n (t) | = 0 
t - * “°° 

for each N, since dim(Zjj) < 00 . 


As in Theorem 2.5, we write 


( 3 * 35 > f n - S 1 B » b n t k’ 


(3.36) E n - R - 


and 


(3.37) S N T n - B n F n . 

From here on, II will be the nonnegative self-adjoint solution to the infinite 
dimensional Riccati algebraic equation (2.14) — when it exists — F will be 
the corresponding feedback operator in (2.18) and S will be the corresponding 
closed-loop operator in (2.21). 
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Theorem 3.9. If II N P N converges strongly to some bounded linear operator II, 
then II is a nonnegative self-adjoint solution to (2.14), F N P N converges in 
norm to F and S N P N converges strongly to S. 

Proof- This follows from Hypotheses 3.1 and 3.2, (3.33) and (3.35) - (3.37), 
and the fact that the control space U has fixed finite dimension. 

Theorem 3.10. Suppose that there exist positive constants M and r, 
independent of N, with r < 1, such that 

(3.38) n N < M, N = 1,2, ... , 
and 

(3.39) |sj| < Mr t , t = 1,2, ... , N = 1,2, ... . 

Then a nonnegative self-adjoint solution II to (2.14) exists, and as 
N ■+■ “, 


(3.40) II P * H strongly. 

If there exists a positive m, independent of N, such that 

(3.41) Q N > m, N = 1,2, ... , 

then (3.38) implies the existence of an r less that one and independent of N 
for which (3.39) holds. 
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Proof* For each N, let II N (- ) satisfy (3.6) with t£ = 0 and II N (0) = where 

I denotes the identity operator on Z^. From (2.32), 

(3.42) |n N - n N (-t) | Has t -► », 
uniformly in N. Now, for z e Z, write 

(3.43) <(n N - n N ,)z,z> z = <(n N - n N (-t))z,z> z + <(n N (-t) - n N ,(-t))z,z> z 

— ^n'^ Z,Z ^Z* 

For £ > 0 choose t > 0 such that I (II. T - II XT (-t))z|„ < E and 
I (II t - n N ,(-t))z| < £. Then, for N and N' large enough, 

| (n^(-t) - n N , (-t))z| z < e. This shows that n^z is a Cauchy sequence in Z 
for each z. Therefore, 11^ converges strongly to a nonnegative self-adjoint 
solution to (2.14). 


An important application of this theorem is when the approximating open- 
loop operators Tjj have an exponential decay rate independent of N, Q is 
coercive and Q N = P n Q| z * t ^ l ^ s case > the zero control gives an upper 

bound, independent of N, on 11^. Such is the case in the example discussed in 
Section 4. 1 and in applications to flexible structures with no rigid-body 
modes and coercive structural damping. 


Theorem 3.11. Suppose that II^P^ converges strongly to II, Q N ? N converges in 
trace norm to Q (hence Q is trace class), and (3.39) holds for positive M and 
r independent of N with r less than one. Then II converges in trace norm 


to n 
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Proof. Since II S^Q^S^II ^ < |s^| H QjsjH i » t ^ e series in (2.26) converges in 
trace norm, uniformly in N. The current result follows then from Lemmas 3.5 
and 3.6. 

Note that Qfj P N converges in trace norm to Q if Q is trace class and 

Q = p n qp nIz n * 

Theorem 3.12 if |JI^| is bounded in N, then a nonnegative self-adjoint 
solution II to (2.14) exists, II^P^ converges weakly to II, and F^Pjj and S^P^ 
converge strongly to F and S, respectively. 

Proof. According to [11 Theorem 6], II^P^ converges weakly to some 
nonnegative self-adjoint bounded II. It follows from (3.33) and Hypotheses 3.1 
and 3.2 that II satisfies (2.14) and that F N and converge as indicated. 

Note that Theorem 3.12 holds if S M P M converges strongly but S P, 

B n N N 

converges only weakly. 

3.3 Implementation of the Approximation Schemes 

In constructing the approximating operators T N , B N ,Q N and G N a standard 
Galerkin approach is often taken; that is, T N = P n T, Bjj = PjjB, 

Qjj = P N^ anc * ^N = P N^* We note however that explicit representations for the 
operators T and B are frequently not available. In particular, this can occur 
when the discrete-time system (2.3) arises from the sampling of an infinite 
dimensional continuous-time system of the form (2.4). In this case it is the 
operators A and 8 which are approximated by a sequence of finite dimensional 
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operators A n and E> n on Z N , from which an approximation to the semigroup 
{ T(s) : s > 0} is obtained as T^(s) = exp (A^s), s > 0. The operators Tjj 
and B n are then T^= T N (t ) and B N = /J T N (s)B N ds, respectively. The strong 
convergence T N P N + T and B^ B is then usually argued using an appropriate 
formulation of the Trotter-Kato theorem, a well known semigroup approximation 
result (see [15] [23]). 

The expressions given by (3.3) - (3.7) are operator equations and 
although they are finite dimensional, they are not appropriate for 
computations. To make use of our approximation framework, we must first 
determine equivalent matrix formulations. Toward this end we assume, without 
loss of generality, that U = R m with the standard basis and inner product and 
let be a basis for Z^. Define the x Gram matrix by 

(3.44) IVy - <+>nV 


For an operator A we denote its matrix representation with respect to the 
bases defined above by [A]. Similarly, for an element z e Z or u e U, we let 
its vector representation be given by [z] or [u] respectively. Standard 
calculations yield 


(3.45) 


* -IT 

[T ] = M [T ] M 
1 N J N 1 N J N 


and 


(3.46) 


* T 

[BJ = [BJ M 


N* 


Defining 
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(3.47) n N (t) = 

(3.4S) Q n = 

and 

(3.49) S„ - M„[G n ] 
we obtain 

(3.50) [uj(t)] = -[F N (t)][z*(t)], t Q < t < t f -1, 

(3.51) [F N (t)] = [R N (t)] _1 [B N ] T n N (t+l)[T N ], 

(3.52) [R N (t)J = [R] + [B N ] T n N (t+l)[B N ], 

(3.53) n N (t) - [T N ] T (n N (t+D - 

n N (t+i)[B N ][R N (t)] _1 [B N ] T n N (t+D)[T N ] + q n , t 0 < t < t f -i, 

(3.54) n N (t f ) = g n . 


Note that since Qjj and are self-adjoint and nonnegative so too are Q 


N 


and G^. Equations (3.50) - (3.54) are therefore in the form of the standard 

ones obtained for the feedback law for a discrete-time linear-quadratic 

K N 

regulator problem in R . Consequently they can be solved using conventional 


techniques. The minimum value of the performance index is given by 
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(3.55) j“ - [4 (t 0 )] T n„(t 0 )t4(t 0 )]. 

Analogously, for the infinite time horizon problem, (3.33), (3.35) and 
(3.36) yield 

(3.56) [u*(t)] = -[F N ][z*(t)], t > t Q , 

(3.57) [F n ] = [R N ] _1 [B N ] T n N [T N ] 

(3.58) [£ n ] = [R] + [B N ] T n N [B N ] 


where 11^ is the solution to the matrix algebraic Riccati equation 


(3.59) n N = [T N ] T (n N - n N [R N ] _1 [B N ] T n N )[T N ] + q n . 


with given by (3.48). The minimum value of the performance index is given 


by 


(3.60) j* - [ z “( Eo )] T n N u"(t 0 )3. 
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4. Examples and Numerical Results 

In this section we describe the application of the general approximation 
framework developed above to a variety of examples. In addition to 
theoretical considerations, in each of the examples below, we discuss some 
numerical results for an inf inite— time horizon problem of the form given in 
Problem (P2). All numerical studies were performed on an IBM Personal 
Computer. The machine we used was equipped with an Intel 8086 math co- 
processor chip and 640K bytes of random access memory (of which less than 384K 
was required). 

Matrix exponentials were computed from eigenvalue-eigenvector 

decompositions obtained using the QR algorithm. The matrix Riccati equations 

(3.59) were solved using a Schur-vector decomposition of the Hamiltonian 

matrix (see [18] [24]). It should be noted that if the eigenvalue pairs of the 

Hamiltonian matrix for a continuous-time linear-quadratic regulator problem 

are asymptotic to ±Y(n) as n -*■ «, then the eigenvalue pairs of the 

Hamiltonian matrix for the corresponding discrete-time problem will be 
. . ±Y(n)x 

asymptotic to e as n + ®. Consequently, for all but very small x, 

conditioning problems arise more quickly than in the continuous-time case when 
the approximating matrix algebraic Riccati equations are solved. 

4. 1. The Heat Equation with Boundary Input 

In this example we consider the scalar parabolic system with boundary 
control given by 

0 3 3 

w(s,x) = a(x) g^- w(s,x), s > 0, x e (0,1) 


(4.1) 
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(A.2) w(s,0) = 0, w(s,l) = v(s) 

(A. 3) w(0,x) = <j>(x) 

with a e H*(0,1), a(x) > 0, x e [0,1], 4> e H^(0,1) = L^CO,!) and v e 1^(0, °°). 

To formulate the discrete-time state equation for this system we let x 
denote the sampling interval and consider only piecewise constant controls v 
given by 

(A. A) v(s) = u(t) s £ [tx, (t+l)x), 

t = 0, 1,2 We choose as our state space Z the Sobolev space 

H^(0,1) with the usual inner product 

(A. 5) <<M> Z = <<M> 0 » jJ* (0)*(0 )d0. 

The state z(t) £ Z is 

(A. 6) z(t) = lim _w(s,*)> t = 1,2, ... 
s+tx 

(A. 7) z ( 0 ) = 4> - 

For t £ {0,1,2...}, we define y(s) £ Z by 


(A. 8) y(s ) =w(s,*) - i|> Q u(t), se (tx, (t+l)x) 
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(4.9) y(tx ) = z(t) — ^qu( t ) > 

where e Z is S iven by ^(x) = x, x e [0,1]. A straight forward 

calculation reveals that y(s) = y(s,») satisfies 

(4.10) y(s) = DaDy(s) + a'u(t), s e (tx, (t+l)x) 

(4.11) y(s)| Q = 0, y(s)| 1 =0, s e (tx,(t+l)x) 

(4.12) y(tx) = z(t) - \J> Q u(t), 

where D denotes the differentiation operator on H^(0,1). 

Let A : dom(A ) c Z + Z be given by 

Dom(A) = H 2 (0,1) n hJ( 0,1) = (<J> e H°(0,1): <|> e H 2 (0,1), 

<(>(0) = <K1) = 0} 

(4. 13) A\f> = DaEty . 

The operator A is densely defined and self-adjoint. It satisfies 

(4.14) <Az,z> z < -<u|z| 2 , z e Dom(A) 

for some u > 0 and has compact resolvent. Also, A is the infinitesimal 
generator of an analytic semigroup of contractions (T(s) : s > 0} on Z which, 
in light of (4.14), satisfies |T(s)| < e , s > 0. It follows therefore, 


that 
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(4.15) y(s) = T(s-tx)y(tx) + J® T T(s-a)a'dau(t) , se [tx, (t+l)x). 
The continuity of y, (4.6), (4.8) and (4.9) imply 

(4.16) z(t) = y(tx) + i|> Q u(t) 
and 

(4.17) z(t+l) = y( (t+1 )x ) + i|> Q u(t), 
and hence that 

(4.18) z(t+l) = y((t+l)x) + ^Qu(t) 

= T(x)(z(t) -* 0 u(t)) + l[ T t+1)T T((t+l)x-c)a'dau(t) + 4> 0 u(t). 
Defining the operators T e L( Z) and B £ L(r\z) by 

(4.19) Tz = T(x)z, z e Z 
and 

(4.20) Bu = [(I - T(t))iJ>q + /q T(a)a'da]u, u e R 1 , 
we obtain 

(4.21) z(t+l) = Tz(t) + Bu ( t ) , t = 0,1,2,.. 
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(4.22) z(0) = 4> . 


We take the performance index to be 

t f~ 1 

(4.23) J(g Q ;0,t f ,<j> ,u) = l {q Q |z(t)|Q + ru(t) 2 } + g Q |z(t f )| 2 

t=0 

with Qq , Sq 0 and r > 0. 

Applying the theory developed in Section 2.1, we have, for the finite 
time interval problem, that the optimal control is given by 

(4.24) u A (t) = -F(t)z*(t), t = 0, 1,2,... ,t f -l, 

where for each t, F(t) is the continuous linear functional on Z given by (2.8) 
- (2.11). It follows that F(t) has a representation f(t,») e H^(0,1) and 
that 

(4.25) F(t)i|> = /J f(t,-)'K0)d0 , 

'P e H°(0, 1 ) , t = 0,1,2,... ,t f -l. 

For the infinite time interval problem (t^ = ", gg = 0), it is 
immediately clear that Hypothesis 2.4 is satisfied. It is also clear that 
(4.14) imples that for each z(0) e H^(0,1), u(t) = 0, t = 0,1,2,... is an 
admissible control, and hence that there exists a unique nonnegative self- 
adjoint solution of the Riccati algebraic equation (2.14). From (2.17) - 
(2.19) we obtain 
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(4.26) u*(t) = -Fz A (t), t = 0,1,2,... 

where F is a continuous linear functional on Z and 

(4.27) Fi|> = /Jf(0)iK0)d9, i|> e H°(0,1), 

with f e H°(0, 1 ). 


We define an approximation scheme using a standard Ritz-Galerkin 
approach. We note that the operator - A is coercive and that it can be written 

& k 

as -A= L L where the operator L and its adjoint L are given by 


(4.28) Up 



+ e hJ(0,1), 


and 


(4.29) 


k 

L = -Da 



¥ e Ao,!), 


respectively. We define the space V = Hq( 0,1) together with the inner 
product 


(4.30) <<M> V = <L<}>,M>> Z , <t>,^ e V. 


We note that V is the energy space associated with the operator - A t 
1 /2 

V = Dom((- A) ' ), and it is the closure of Dom(A) with respect to the energy 
norm |»| , which satisfies 



-36- 


(4.31) |*|J = <L4» ,14>> z = <— A <f> ,<t> > z , 


<j> e Dom(A). 


For each N = 2,3,... let A^ denote the uniform partition of the interval 
[0*1] given by {0, , ...,^p- , 1}. Let denote the usual linear 

B-splines on [0,1] corresponding to the partition A^ and which satisfy 
e^(0) = e ^(l) = 0» j = 1,2,... N-l. The ejj are given by 


/ N(0 - 

(4.32) ej (0) = N(-£J±Ii - 0) 
l 0 


9 E [J F- • N 1 

0 s £ • J S L I 


elsewhere 


j = 1,2,..., N-l and are depicted in the figure below. 



Figure 4.1 

i N-l 

Let Z^ = span{e^}j_j . We note that Zjj c V and dim Z N = N-l, N = 2,3,... . 
Define P^:Z + to be the orthogonal projection of Z onto Zjj with respect to 
the < , >2 inner product and : V ->■ Z^ to be the orthogonal projection of V 
onto Zjg with respect to the < , >y inner product. 

We define the operator : Z N Z^ as the inverse of the operator given 
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by 


(4.33) 




J N 


The invertibility of A of course follows from the coerclvity of -A (see 
(4.14)). On the other hand, a straight forward calculation yields 


(4.34) < A n z n , z n > v l z N l z * Z N G Z N* 


Consequently the operator A ^ is invertible and the operator A^ is well 
defined. It is also self-adjoint. Indeed 


(4.35) <A N z N’ Y N > Z " ^N’^V ’ Z N’ Y N G Z N* 

From (4.35) it also follows that 

(4.36) < I Z N ^Z ’ Z N G Z N* 

It can be concluded therefore, that A N is the infinitesimal generator of a 
semigroup of bounded linear operators on Z N , { "^(s) : s > 0} with 
T n (s) = exp(A N s), s > 0 and satisfying 

(4.37) |T n (s)| z < e^ s , s > 0. 

Elementary properties of spline functions (see [27]) imply I 

as N + strongly on V. Furthermore, A * compact and 
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(4.38) 




|(P N - I) A -1 z | y 


imply that P N A + A as N + ® in the uniform operator topology on L(Z). 
We have therefore, that 


<4 ' 39) l A N lp N ~ A_1 | Z * ° 

as N + 

From (4.37) and (4.39) we conclude (see [4], [15]) 

(4.40) T N ( s )P N z T(s)z 


and 


(4.41) T n (s)P^z T (s)z 

as N 00 for each z e Z, uniformly in s for s in compact subintervals. 

Wlth T n ' V T) - % ' W G N ■ «0 P N “<1 

(4.42) B N - (I - T|) )P n (. + ! r 0 T 1) (a)P N a'do, 

(4.40), (4.41) and elementary approximation properties of spline functions 
guarantee that Hypotheses 3.1 and 3.2 hold and hence that the convergence 
results for the finite time interval problem given in Theorem 3.3 apply. 

For the infinite time interval Problem, (4.37) implies that Hypothesis 
3.8 is satisfied. Moreover, if qg > 0 the conditions given in the statement 
of Theorem 3.10 are satisfied (see the Remark following the proof of Theorem 3.10) 
and consequently the convergence results for the infinite time horizon problem 
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given in Theorem 3.9 hold. 

Define the ^ vector valued function E M on [0,1] by 


(4.43) E n ( 0) T = (e*(0), ej(0),..., e* _1 (0)). 


0 e [0,1] and the (N - 1) * (N - 1) Gram matrix associated with the basis 

, r ii N_1 

elements ie. T l . , 

N j = l 


(4.44) 


M N = <E 


N* 


T 

Vz 


A straight forward calulation yields 


(4.45) 



Let the (N-l)x(N-l) matrix H N be given by 


(4.46) H n - <E„, . 

Then 


iV u 



A 0 e N 


> 


z 



e j > 

N V 


(4.47) 
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<Le N , Le^> z 


' ‘ <aD 4- De N > Z • 


In the case a(x) - a, x £ [0,1], a constant, we have 



The matrix representation [ A N ] for the operator A N is given 


(4.49) 


M N H N* 


while for the operators T N , Q N and G N we have 


[T N ] = ex PdA N JO> 

[ V = q 0 I N 

and 


(4.50) 

(4.51) 


(4.52) 


[ V g 0 I N 
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where I N denotes the (N-l)x (N-l) identity matrix. If we define 
♦0N> a N E RN " 1 b? 

(4.53) * 0[) - <E n , Vz 


and 

(4.54) a^ = <E N ,a’> z , 
j = 1,2,..., N-l we obtain 

(4.55) [B n J = (I N - [T n ])^ on + J] exp([A N ]a)M" 1 a^da 

= ^ON + a N * 

When the finite dimensional approximating gain matrices [F N (t)], t = 
0,l,2,...t f -l for the finite time interval problem have been computed using 
(3.51)-(3.54), an approximation f N (t,«) to the feedback kernal f(t,») can be 
obtained from 

(4.56) f N (t ’ 0) = V 9)TM N 1[F N (t)]T ’ 
t = 0, 1 ,2, . . . t^-1 , 9 e [0,1]. We have 

f N (t,«) * f(t,») 


(4.57) 
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in H®(0,1) as N + “ for each t = 0, 1 ,2, . . . t^-1. 

Similarly, for the infinite time interval problem, an approximation f N to 
f is given by 

(4.58) ( N (e) = E s (e)\ 1 [F N ] T , 

0 e [0,1] where the matrix [Fjj] is computed from (3.57) - (3.59). We have 

(4.59) f N f 
in H° ( 0 , 1 ) as N + » . 

We demonstrate the feasibility of our schemes on an infinite time 
interval problem of the form discussed above. Taking qg = 1.0, r = 1.0, 
a(x) = a = 1.0, x £ [0,1] and x = .01 we obtained the approximating feedback 
kernals given in Table 4.2 and shown in Figure 4.3 below. 



3 - 


-4 



f 3 (6 > 

f 5 ( 0 ) 

f 9 ( 0 ) 

fu «0 

fl3<0) 

0.00 

.0000 

.0000 

.0000 

.0000 

.0000 

0.05 

.0250 

.0228 

.0224 

.0222 

.0221 

0.10 

.0500 

.0456 

.0449 

.0445 

.0443 

0.15 

.0750 

.0684 

.0673 

.0668 

.0664 

0.20 

.1000 

.0911 

.0898 

.0891 

.0887 

0.25 

.1250 

.1147 

.1122 

.1115 

.1109 

0.30 

.1500 

.1382 

.1346 

.1338 

.1331 

0.35 

.1727 

.1617 

.1571 

.1561 

.1554 

0.40 

.1908 

.1853 

.1798 

.1787 

.1778 

0.45 

.2088 

.2129 

.2028 

.2013 

.2005 

0.50 

.2269 

.2405 

.2286 

.2261 

.2245 

0.55 

.2450 

.2681 

.2543 

.2513 

.2496 

0.60 

.2630 

.2958 

.2819 

.2793 

.2769 

0.65 

.2811 

.2956 

.3097 

.3056 

.3035 

0.70 

.2584 

.2954 

.3246 

.3275 

.3275 

0.75 

.2153 

.2953 

.3331 

.3352 

.3385 

0.80 

.1723 

.2951 

.3169 

.3260 

.3263 

0.85 

.1292 

.2213 

.2697 

.2163 

.2170 

0.90 

.0861 

.1475 

.2097 

.2163 

.2170 

0.95 

.0431 

.0738 

.1049 

.1125 

.1175 

1.00 

.0000 

.0000 

.0000 

.0000 

.0000 


Table 4.2 
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4.2. Hereditary or Time Delay Systems 

In this example we consider linear hereditary control systems of the form 

(4.60) x(s) = Lx + B~v(s) s > 0 

s U 

where x(s) £ R n , x g £ L2((-r,0);R n ) for some r > 0, v £ L2( (0, s^ ) ;R m ) for all 
s < “ and Bq e L(R m ,R n ). The function x g represents the history on the 
interval [s-r,s]; that is x g (0) = x(s+6), 0 £ [-r,0]. The operator L is 
assumed to be of the form 

v 

(4.61) L<j> = l A <|> (-r ) + /_ A(0 )<j) (0 )d0 

i=0 

where A.^ £ L(R n ), i = 0,l,2,...,v, A(» ) £ L2((-r,0); f(R^)) and 
0 = Tq < r^ < ••• < r y = r. The initial data is given by 

(4.62) x(0) = n , x Q = <f) , 

where n e R n and <|> £ L 2 ((-r,0); R°). 

Once again, we formulate the discrete-time control problem by letting x 
denote the sampling interval and considering piecewise constant control inputs 
of the form 

(4.63) v(s) = u(t) s £ [tx, (t+l)x), 

t = 0,1,2,..., where u(t) £ R m for each t. The state space is 
Z = R n x L 2 ((-r,0);R n ) 


(4.64) 
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with the inner product 

(4.65) <(n ,4>), (5 ,tP>„ = <n,C> + <<M> T . 

z R n l 2 

The state of the discrete-time control system is 

(4.66) z(t) = (x(tx ;n ,4> ,v) , x (n ,<|) , v) ) , t = 0,1,2,..., 

where x(s;n,<t>,v) denotes the solution at time s to the system (4.60) - (4.62) 
and x (n,<f> ,v) its history on [s-r,s]. Then 

(4.67) z(t+l) = T(t) z(t) + Jq T(o) Bdau(t), t = 0,1,2,..., 

(4.68) z(0) = (n ,<(> ) 

where (T(s) : s > 0} is the Cq semigroup of bounded linear operators on Z 
with infinitesimal generator A: Dom(A) c z -*■ Z given by 

Dom ( A ) = {(£,<!>) e Z : e H 1 ((-r,0); R n ), ? = ij>(0)} 

(4.69) A0K0), 4 .) = (Lt|», D>J»), 

and the operator B : R™ ^ Z is defined by 

(4.70) Bu = (Bqu, 0), u g R m . 

Letting T = T(x) and B = /g T(o)Bdo, we obtain 
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(4.71) 


z(t+l) = Tz(t) + Bu(t) , t = 0,1,2,... 
z(0) = (n,<j)). 


Let Qq, Gq e L(R n ) be symmetric and nonnegative and let R e L (R m ) be symmetric 
and positive definite. We consider the discrete-time linear-quadratic 
regulator problem with state given by (4.71) and performance index 


(4.72) J(G;0,t f ,n,<j>,u) 


1 {<Qz(t),z(t)> + 

t=0 L 


<Ru(t) ,u(t)> } + <Gz(t ),z(t^)> 

„m f f Z 


where the operators Q : Z Z and G : Z -»■ Z are given by 


(4.73) Q(5,¥) = (Q q 5,0) 


and 

(4.74) G(e,1>) = (G q 5,0) 

respectively. For the infinite time problem we of course have 
t = * and G q = 0. 

For the finite time problem with the operators, T,B,Q and G as defined 
above, we apply the theory developed in Section 2.1 and conclude that the 
optimal control is given by 


t 0,1,2,..., t 1 , 


(4.75) 


u*(t) = -F(t)z jl[ (t), 
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where the operators F(t) £ L(Z,R n ) are given by (2.8) - (2.11). The operator 
F(t) can be represented by a matrix of operators, [F®(t), F*(t)] with 
F°(t) e l(R n ,R m ) and F^t) e L (L 2 ( (-r,0); R n ); R m ). It follows from (4.70) 
therefore that 

(4.76) u A (t) = -f°(t)x A (tx) - /° r f 1 (t,0)(x^) tx (0)d0, t = 0,1,2,... ,t f -l, 

where f^(t) is an m x n matrix and f*(t,*) is a square Integrable m x n 
matrix valued function on (-r,0). 

For the infinite time problem we assume that our original hereditary 
system and Qq are such that there exists an admissible control for each 
z(0) = (n,<{>) £ Z and that Hypothesis 2.4 is satisfied. Then Theorems 2.3 and 
2.5 imply that there exists a unique nonnegative self-adjoint solution II to 
the Riccati algebraic equation (2.14) with the optimal control u* given by 

(4.77) u*(t) = ~Fz*(t), t = 0,1,2 

where F £ L(Z,R m ) is given by (2.18) - (2.19). The feedback gain F can be 
represented by a matrix of operators [F^,F*] with F^ £ L(R n ,R m ) and 
F^ £ L(l> 2 ( (-r ,0) ;R n ) ,R m ). We have 

(4.78) u*(t) = -f°x*(tT) " /° r f 1 (0)(x*) tx (0)d0, t = 0,1,2,..., 

where f^ is an m x n matrix and f* is a square integrable 
m x n matrix valued function on (-r,0). 

Numerical methods for the approximate solution of the continuous time 
linear quadratic regulator problem for hereditary systems in closed-loop form 
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have been studied extensively (see [7 ] , [ 16 ] , [ 19 ] , [25 ] , [26 ] ). Most closely 
related to the approximation framework which we have developed here for the 
discrete-time problem are the treatments for the continuous-time problem in 
[2] , [ 11 ] , [ 14] and [17]. The first approximation scheme applied to the 
continuous-time linear-quadratic control problem for hereditary systems was 
the AVE scheme used in [1 ] , [ 11 ] , [25] and [26] which approximates the history 
by piecewise constant functions and derivatives with finite differences. In 
[1], Banks and Burns cast AVE in its modern semigroup approximation form, 
proved strong convergence of the approximating open-loop semigroups and used 
the scheme to compute the open-loop optimal control. Gibson [11] used the 
convergence results in [1] to prove strong convergence of the solutions to the 
approximating continuous-time Riccati equations and uniform norm convergence 
of the approximating feedback control laws. However, the rate of convergence 
as it was observed in numerical studies was relatively slow. Spline based 
schemes for continuous- time systems were developed for the finite time 
interval problem in [17] and for the infinite time interval problem in [2]. 
With regard to the minimization of the performance index these schemes 
represented an improvement over the AVE scheme. However, they appeared to 
yield only weak convergence of the solutions to the approximating Riccati 
equations and strong operator convergence of the approximating feedback gains. 

Recently, a new spline-based state approximation for hereditary systems 
has been proposed in [14]. When applied to the continuous time control 
problem, this new method performs at the level of the earlier spline-based 
schemes and yields strong convergence of the approximating feedback gains. We 
have chosen this scheme to describe and implement here for the discrete-time 


problem. 
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To simplify the presentation, we consider systems of the form (4.60) 
having only a single discrete delay (v = 1) and no distributed delay term 
(A(0) = 0). A detailed description of the scheme in complete generality can 
be found in [14], 

For each N = 1,2,... define 

(4* 79 > e N = (I n ,0) and ej = (O.eJ ), j = 1,2,. ..N+l, 

where I n denotes the n * n identity matrix and the e^ are the 
n x n matrix valued piecewise linear elements on [-r,0) given by 


7 ( 0 + £)i 

r N 


n 


eJC 0 ) - 


0 e II . o) 

elsewhere 


(4.80) 


( 


7 ( 0 + J £ )i 


n 


ej( 0 ) ={ 


-N 

r 


? (9 + (j-2) f)I r 


\0 


e e H f , -(3-1) £) 

0 e [-(J-l)f , ~(j-2)f) 

elsewhere 


j = 2,3,4,.. .N 
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‘f (9 + (N-l) X 


eT<9> - 


0 e [-r,-(n-l) — ) 


elsewhere . 


Let 


N+l 


(4.81) Z. T = {z e Z : z = 7 a . e^ T , a. e R n }. 

N jto J N J 


Note that dim Z^ = n(N +2). We shall refer to the collection { as a 


N+l 


'N j=0 


"basis" for Z„ and a vector a e x r“ as being a "coordinate vector" for an 


J N 


element in Zjj. Defining 


j=0 


f, QON /0 "1 "N+l. 

(4.82) E n - <e N , e„,— e N ) 


we have 

. N+l 

(4.83) Z = {z e Z : z = E„a, ct e x R n } . 

” " J-O 

Let Mjj denote the Gram matrix corresponding to the basis 
straight forward computation yields 


(4.84) M n = diag (I n , ^ ®I n ) 


where the (N+l) x (N+l) matrix m^ is given by 
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(4.85) 



and ® denotes the Kronecker product. 

Let Pjj denote the orthogonal projection of Z on to Z N with respect to the 
inner product (4.65). It follows that 


(4.86) P n (5,i|0 = (5,p n 10 

where p N is the orthogonal projection of L 2 ((-r,0) ;R n ) on to spante^}*?** with 
respect to the usual L 2 inner product. We have 

(4.87) P n (5,i|0 = £^(5, <10 

N +1 

where a (C,')') e x R n is given by 
j=0 

(4.88) «„«,*) -M' 1 col «,♦*, 


with 


K = f-r 4 (B)i,(6)de - 


(4.89) 



-53- 


We shall set = exp( A^x) where A^ : -► is an appropriately 

defined finite dimensional approximation to the operator A given by (4.69). 
Noting that Z N t- Dom(A), we motivate the definition of A N by first formally 
extending the operator A to an operator defined on Z^. 

F ° r Z N = ( W £ Z N d6fine 

(4.90) K\ - (A^+A^-r), D + + N + S<? N - lii»_t H (0 ))) 

0+O 

where 6 is the Dirac delta impulse concentrated at zero and D% denotes the 
right hand derivative of For each N = 1,2... let A N : Z N •>• Z N be given by 

(4.91) A nI „ - (VAt,!-!), p/V + «„«„ - lta.ye)) 


where 


(4.92) 


6 n 


;t 

e n y n 


with 

(4.93) Y = M 1 col (0, lim_ e*(0), ••• , lim e!! +1 (0)). 

0+O 0-K) - 

To compute [A N J, the matrix representation for the operator A^, we let 
[z N ] denote the coordinate vector representation for an element z N e Z N « 
Then from 


(4.94) [A n z n ] = [A n ][z n ] 
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and 


(4.95> [A n z n ] 


H ti a tJ (A N Z tI ) m n <e n ,a n z n > z m n <e n ,a n e b > z^ z n' 


we obtain 


(4.96) [A ] = m" 1 H 
N N N 


where 

(4.97) H n = « n ,A^> z . 


Using the definitions of A^ and a straight forward calculation (see 
[14]) yields 


(4.98) 



where the (N+l) x (N+l) matrix hjj is given by 
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(4.99) 



The matrix representation for the operator T N = exp (A^t ) can then be 
computed from 


(4.100) [T n ] = exp([A N ]T). 


We define the operators : R™ * Z. T , Q„ : Z' T Z XT and G„ : Z„ ■»■ Z„ by 

N N NNN NNN 


(4.101) B N = P N 15 

(4.102) Q n = p n Q 

and 

(4.103) G n = P n G 


from which we obtain 
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(4.104) 


[8 n ] = col (B 


0 ’ 


0,... ,0) £ R 


n(N+2)xm 


(4.105) 



and 


(4.106) [G n ] 



0 


n(N+2)x n (N+2) 

G K 


o 


Finally, defining 

(4. id 7 ) B n = Jj exp (A^s) B N ds 


we have 

(4.108) [B n ] = /q exp ( IA N ls) [ B^lds. 


Once the matrix representations for the approximating feedback gains have 
been computed, [F N ( t ) ] , t = 0,1,2,... t f -l from (3.51) - (3.54) for the finite 
time interval problem and [F N ] from (3.57) - (3.59) (assuming, for the moment, 
that solutions to (3.34) exist) for the infinite time interval problem. 
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approximations for f^, f*(t,*)» t = 01,2,. ..t f -l and can be 

computed from 


(4.109) ((f°(t)) T , (fJ(t,.)) T ) -^ 1[ Vt)l T , 


t = 0,1,2,. 





and 


(4.110) ((f°) T , (f°(*)) T ) ^^[F/ 


respectively. 

For the approximation scheme defined above, it is shown in [14] that 
P N * I strongly on Z. Using a Trotter-Kato like result it is also shown that 

(4.111) exp (A^s)P^ T (s) 


and 

(4.112) exp(A N s)P N -»■ T (s) 

strongly on Z and uniformly in s for s in compact intervals. Hypothesis 3.1 
is a simple consequence of these results. The present scheme, therefore, 
satisfies all of the hypotheses of Theorem 3.3 and we may conclude that the 
convergence results for the finite time interval problem given in the 
statement of the theorem hold. In particular, we have 

(4.113) f°(t) ■* f°(t) 
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. „mxn , 
in R and 

(4.114) f*(t,«) * f^t,*) 

in L 2 ((-r, 0 );R nf<n ) for each t = 0, 1 ,2, . . . t^-1. 

With the operators Q and G given by (4.73) and (4.74) and the operators 

Qjj and defined as in (4. 102) and (4. 103) it is clear that the hypotheses 

given in the statement of Theorem 3.7 are satisfied. We have therefore that 

*•£ 

for the present example the operators {II(t)} t _Q are trace class and 

(4.115) lim in (t)P - n ( t ) II = 0, t = 0,1,2,. ..,t_. 

N-n» N IN i t 

For the infinite time problem and the approximation scheme discussed 
here, the situation with regard to convergence is much the same as it is for 
the continuous time problem (see [14]). We are unable to demonstrate the 
existence of an M and an r < 1 for which (3.38) and (3.39) hold. In fact, our 
numerical studies point to the conclusion that condition (3.39) is violated by 
the present scheme. We observe the existence of a sequence of closed-loop 
eigenvalues of the approximating discrete-time control problems (p2jj) which 
tend toward the unit circle as N -*■ °°. 

On the other hand, upon solving the approximating problems it is also 
apparent that |ll^| remains bounded in N. Consequently we may apply Theorem 
3.12 to conclude that a solution II to (2.14) does in fact exist, 

HjgPjg + n weakly and -> F strongly as N + » . We have therefore that 


f 


0 

N 


f 


0 


(4.116) 
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, „mxn 
in R and 

(4.117) f N ^ fl 

weakly in L^C (-rjCOjR 11 ^ 11 ) as N + 00 . 

We applied the scheme to the infinite-time problem with state 


(4.118) y(s) + y(s-l) = v(s ). 


Transforming (4.118) to an equivalent first order system we obtain a system of 
the form given in (4.60), (4.61) with n = 2, r = 1, m=l. 


’o 

il To 

oi r 

ol 

' y(s) 

0 

°J* A i |_ _1 

O 1 
Cd 

o 

II 

jl, A = 0 and x(s) = 

_y(s). 


Taking Qq - |^q the performance index takes the form 


(4.119) J(0;0,oo,y(0),y(0),y 0 ,y 0 ,u) 


- I 

t=0 


y( tx )' 


(tx) 2 + 


Ru(t)‘ 


where t is the length of the sampling interval. The optimal feedback control 
is given by 


u*(t) 

(4.120) 


-[f°] iy (tT) -[f°] 2 y(tT) 

- /° 1 {[f 1 (0)] 1 y(tT+0) + [f 1 (0)] 2 y(tx+6)}d0 


where [f°] 1 and [f 1 (0)]_ L , i = 1,2 are the i th components of the 1 x 2 
matrices f^ and f*(0). 
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Since by taking the initial conditions 
(4.121) y(0) = 0, y(0) = 0, y o (0) = 0, -1 < 0 < 0 

we have y(s) =0, s > 0 regardless of how yg(0), -1 < 0 < 0 is chosen, it 

follows that [ f ^ (0 ) ] 2 = 0, -1 < 0 < 0. Indeed, the optimal control 

corresponding to the initial conditions (4.121) with y^ arbitrary is u(t) = 

0, t = 0,1,2 Furthermore, the nature of the approximation scheme is 

such that we must have [f*(0)] o = 0, -1 < 0 < 0, N > 1. 

N 2 

Setting x = .01 we obtained the results given in Tables 4.4 and 4.5 and 
Figure 4.8 below when R = .05. With R = 1.0, the results given in Tables 4.6 
and 4.7 and Figure 4.9 were obtained. As the cost of control increases the 
effect that the optimal control for the infinite dimensional problem has on 
higher modes decreases. Consequently, the finite dimensional approximations 
are more effective and convergence is more rapid. 
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N 

2 

4 

8 

10 


4.5483 

4.5452 

4.5451 

4.5451 

*- f N ^2 

5.2954 

5.2948 

5.2948 

5.2948 


Table 4.4 


0 

warn 

»!«>>]! 


IfJoCfOli 

.00 

.1277 

.0942 

.0878 

.0872 

-.05 

.0948 

.0690 

.0659 

.0658 

-.10 

.0619 

.0437 

.0439 

.0445 

-.15 

.0289 

.0185 

.0163 

.0154 

-.20 

-.0040 

-.0068 

-.0171 

-.0137 

-.25 

-.0369 

-.0321 

-.0506 

-.0556 

-.30 

-.0698 

-.1046 

-.1054 

-.0976 

-.35 

-. 1027 

-.1772 

-. 1603 

-.1650 

-.40 

-.1357 

-.2497 

-.2375 

-.2324 

-.45 

-.1686 

-.3223 

-.3371 

-.3361 

-.50 

-.2015 

-.3949 

-.4367 

-.4399 

-.55 

-.5891 

-.6156 

-.6090 

-.6054 

-.60 

-.9767 

-.8363 

-.7813 

-.7709 

-.65 

- 1.3643 

- 1.0570 

- 1.0200 

- 1.0283 

-.70 

- 1.7519 

- 1.2777 

- 1.3251 

- 1.2858 

-.75 

- 2.1395 

- 1.4984 

- 1.6301 

- 1.6899 

-.80 

- 2.5271 

- 2.1836 

- 2.1603 

- 2.0941 

-.85 

- 2.9147 

- 2.8688 

- 2.6905 

- 2.7232 

-. 90 

- 3.3023 

- 3.5540 

- 3.4158 

- 3.3523 

-.95 

- 3.6899 

- 4.2391 

- 4.3361 

- 4.3269 

1 

• 

o 

o 

- 4.0775 

- 4.9243 

- 5.2565 

- 5.3017 


Table 4.5 


















- 62 - 


N 

2 

4 

8 

10 

l£ S ] i 

1.4050 

1.4054 

1.4054 

1.4054 

“£>2 

1.9477 

1.9479 

1.9479 

1.9479 


Table 4.6 


0 

[f 2( 0 )] 1 . 


■ssatgak si 

i£ io (e)i i 

.00 

-.2813 

-.2831 

-.2835 

-.2835 

-.05 

-.3433 

-.3402 

-.3383 

-.3379 

-.10 

-.4052 

-.3972 

-.3931 

-.3923 

-.15 

-.4052 

-.4543 

-.4522 

-.4534 

-.20 

-.5292 

-.5114 

-.5156 

-.5146 

-.25 

-.5911 

-.5684 

-.5790 

-.5798 

-.30 

-.6531 

-.6439 

-.6478 

-.6450 

-.35 

-.7151 

-.7195 

-.7165 

-.7179 

-.40 

-.7770 

-.7950 

-.7902 

-.7907 

-.45 

-.8390 

-.8705 

-.8689 

-.8685 

-.50 

-.9010 

-.9460 

-.9476 

-.9462 

-.55 

- 1.0044 

- 1.0347 

- 1.0328 

- 1.0324 

-.60 

- 1.1078 

- 1.1233 

- 1.1181 

- 1.1185 

-.65 

- 1.2112 

- 1.2120 

- 1.2089 

- 1.2103 

i 

• 

O 

- 1.3146 

- 1.3007 

- 1.3053 

- 1.3021 

-.75 

- 1.4180 

- 1.3894 

- 1.4016 

- 1.4030 

-.80 

- 1.5214 

- 1.5013 

- 1.5057 

- 1.5040 

-.85 

- 1.6248 

- 1.6132 

- 1.6098 

- 1.6112 

-. 90 

- 1.7282 

- 1.7251 

- 1.7198 

- 1.7185 

-.95 

- 1.8316 

- 1.8370 

- 1.8357 

- 1.8352 

- 1.00 

- 1.9350 

- 1.9489 

- 1.9516 

- 1.9519 


Table 4.7 
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4.3 Control of a Flexible Structure 

We consider an Euler-Bernoulli beam cantilevered to a rigid hub which is 
free to rotate about its fixed center, point 0. Also, a point mass m^ is 
attached to the other end of the beam. The control is a torque u applied to 
the hub, and all motion is in the plane. See Figure 4.10 and Table 4.11. 



u(t) 


Figure 4. 10 
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r = hub radius 
H = beam length 

Iq = hub moment of inertia about axis 

perpendicular to page through 0 
m^ = beam mass per unit length 
m^ = tip mass 

El = product of elastic modulus and 

second moment of cross section for beam 
fundamental frequency of undamped structure 


10 in 
100 in 


100 slug in^ 

.01 slug/in 
1 slug 

13,333 slg in^/sec^ 
.9672 rad/sec 


Table 4.11 


The angle 0 represents the rotation of the hub (the rigid-body mode), 
w(t,n) is the elastic deflection of the beam from the rigid-body position, 
and wj(t) is the displacement of m^ from the rigid-body position. For 
technical reasons, we do not yet impose the condition Wj(t) = w(t,£). 

The control problem is to stabilize rigid-body motions and linear (small) 
transverse elastic vibrations about the state 0=0 and w = 0 . Our linear 
model assumes not only that the elastic deflection of the beam is linear but 
also that the axial inertial force produced by the rigid-body angular velocity 
has negligible effect on the bending stiffness of the beam. The rigid-body 
angle need not be small. 

For this example, it is a straight forward exercise to derive the coupled 
ordinary and partial differential equtions of motion in 0 , w and w^. However, 
rather than writing these equations explicitly, it is easier and more useful 



for our purposes to derive an abstract second order evolution equation for the 
structure. To do this, we define the generalized displacement vector 


(4.122) 


x = (9,w,Wj) £ H = R* x 1^(0, £) x R*. 


The kinetic energy in the system is then 


(4.123) Kinetic Energy =1/2 <MqX,x>^ 


where Mq is the unique bounded self-adjoint linear mass operator Mq on H such 
that 


(4.124) <MqX , x> = I Q 00 + m b <w + ^ Q 0,w + ^ O 0> L + m l^ w i + + ^q^ 9 ^’ 


where i|>q £ 1^(0, £) is given by 'I'q('I) = r + n. It is easy to show that Mq is 
also coercive. The elastic strain energy is 


(4.125) Strain Energy = 1/2 a(x,x) 
with 

(4.126) a(x,x) = El <D^w, D^w> T . 


We make a(*, # ) into an inner product by setting 
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and define the strain-energy space 

(4.128) V = (x = (0 ,<M (£)): <(> e H 2 (0,i), <J> (0) = Dc{> (0) = 0} . 


The last term in (4.127) is necessary for the V-inner product because there is 
no strain energy associated with the rotation of the hub. 


We define the stiffness operator Aq by 
(4.129) Dom(A Q ) = {x = (0 ,♦,♦(£)) e V: $ e H 4 (0,£), D 2 <j>(£) = 0} 


and 


(4.130) 



0 0 0 

4 

0 El D 0 

0 -El D 3 0 


This operator is self-adjoint with compact resolvent and all positive 
eigenvalues except the one zero eigenvalue corresponding to the rigid-body 
mode. Note that V is the domain of the square root of Aq. 

With these mass and stiffness operators, we can write the equations of 
motion as 

(4.131) Mqx(s) + CqAqx(s) + Aqx(s) = B Q u(s), s > 0, 

where Cg is a positive constant and the term CqAqX represents viscoelastic 
damping in the beam. The input operator is 
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(4. 132) B q = (1,0,0). 

A A A A 

Letting Z = V x H with inner product <(v,h) , (v,h)> z = <v,v> v + <K Q h,h> H , 
the first order form of this system is given by 

(4.133) z(s) = Az(s) + Bu(s), s > 0, 


where z = (x,x) £ Z and A is the unique extension of the operator 


(4.134) 


O 

A = 


L M 0 A 0 " M 0 A 0 J 


o 

Dom(A) = Doiii(Aq) x Dom(Ag), 


that generates a CQ-semigroup on the space Z. Of course, 8 is 


(4.135) 


8 = 


l m o 1b o 


See [10] and [12], The hub-beam-tip mass structure here is discussed in more 
detail in [12], along with the continuous-time problem. 

The discrete-time control system for sampling interval t is 


(4.136) z(t+l) = Tz(t) + Bu(t), t= 0,1,2,..., 


where 

(4.137) T = T (x ) , B = /q T(s)B ds 


and (T(s): s > 0} is the semigroup generated by the A in (4.133) 
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As in the previous examples, we will solve a discrete-time optimal 
control problem on the infinite interval. In the performance index, we take 
the state weighting operator Q to be the identity on Z. This means that 
<Qz,z> z is twice the total energy in the structure plus the square of the 
rigid-body rotation. Since there is one input, the control weighting R is a 
scalar. The optimal control has the feedback form 

(4.138) u*(t) = " <f,x(t)> v - <M Q g,x(t)> H 

where x(t) has the form (4.122) and 

(4.139) f = (f\f 2 ,f 3 ) e V 

(4.140) g = (g^g^g 3 ) e H . 


Our approximation of the structure is based on a finite element 

approximation of the beam which uses Hermite cubic splines as basis functions 

J 

([27]). We define the sequence of spaces V N = span {e^}j_^ with 


(4.141) 


ej = (1,0,0), 


(4.142) 


e N «>.♦»»■ ^ 2,3,... J^, 


where the c^'s are the cubic splines. Each Vjj is a subspace of V, and our 
approximation scheme is a Ritz-Galerkin approximation obtained by projecting 
(4.131) onto Vjj. See [12] for details. Writing 
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(4.143) 


x N (s) 


l lV s)I j e K> 

J=1 J 


we have 


(4.144) 


Mjj [^(s)] + CqK^ [x n (s)] + K n [x n (s)] BqjjU(s) 


to solve for the vector [x N (s)] of time-dependent coefficients [x N (s)]j. The 
mass matrix and the stiffness matrix are given by 


(4.145) 


[MJ 


N ij 


= <M e e J > 

^ 0 e N’ TH’ 


* K N J ij <e N* e N > V 


and the input matrix is 


(4.146) B Qn =[100 ... 0] T . 

With z^ = (x^,x N ) e V N x V (4.144) is the matrix representation of an 
evolution equation 


(4. i47) z N (s) =A n z n (s) + B n u(s) 

where A N and B N approximate A and 8 . It is shown in [12] that, as N 

increases, the semigroup {T^(s): s > 0} generated by Ajj converges strongly to 

the semigroup {T (s) : s > 0} and that the adjoint semigroup (T*(s): s > 0} 

converges strongly as well. Since 8„ is the Z-projection of 8 onto V x V , 

N N 

it converges strongly to B . 

For the approximating discrete-time control systems, we replace z(t), T, 
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B, T( m ) and B in (4.136) and (4.137) with z N ,T N ,B N , T (• ) and B N , 
respectively. For each N, the solution to the infinite-time optimal control 
problem is based on the N c Riccati operator equation (3.33). As in the 

a 

previous examples, we solve the Riccati matrix equation (3.59) for II N , which 
is related to [11^] (the matrix representation of the operator II N ) as in 
Section 3.3, except here we have 


(4.148) 


fl N W N [n N ] ’ 


where 


(4.149) 


W.. 




and is the stiffness matrix with 1 added to the first element. Since Q = 

I in the infinite dimensional problem, Q M is the identity on V x V and it 

w N N 

A 

follows from (3.48) that the matrix for (3.59) is W^. 

The optimal feedback control for the N 1 "* 1 problem is given by (3.56) with 
the matrices in (3.57) and (3.58), and it has the equivalent representation 


(4.150) 


N, . 
u *(t) = 




0 S N’ N 


H 


where 


(4.151) 


f N ^ f N ,f N ,f N^ E V ’ 


(4.152) 


, 1 2 3. _ 

g N ^ g N ,g N’ g N^ £ H 



-72- 


as in (4.139) and (4.140). From (3.56), (4.143) and (4.150), it follows that 



the control weighting R = 1, Tables 4.12-4.15 give the values of the 

corresponding scalar and functional gains, f^, g*, i = 1,2,3 for various 

values of N. The values of the functional gains D and g^ along the length 

of the beam also are plotted in Figures 4.16 and 4.17. We tabulated and 
2 2 

plotted D f because this is what appears in the V inner product in (4.150) 

2 

and also to show the H convergence. We note that the form of the V inner 

3 

product given by (4.127) is such that f^ does not appear in the feedback law. 




Table 4.12 


N 

3 

4 

5 

7 

f N 

.1750 

.1769 

.1774 

.1777 

3 

% 

-18.1231 

-18.3385 

-18.3902 

-18.4158 


Table 4.13 























0.0 


6.0266 


7.0718 


7.6281 


8.1441 


5.0 

4.8991 

5.4251 

5.6142 

10.0 

3.7717 

3.7784 

3.6003 

15.0 

2.6442 

2.1317 

1.5865 

20.0 

1.5167 

.4850 

.2048 

25.0 

.3893 

.1353 

.4057 

30.0 

-.7382 

.4896 

. 6066 

35.0 

1.0908 

.8438 

.8075 

40.0 

1.2767 

1.1980 

1.3790 

45.0 

1.4626 

1.5523 

1.5713 

50.0 

1.6484 

2.0984 

1.7636 

55.0 

1.8343 

1.9858 

1.9559 

60.0 

2.0202 

1.8731 

2.0543 

65.0 

2.2061 

1.7605 

1.8179 

70.0 

1.6333 

1.6478 

1.5815 

75.0 

1.3545 

1.3288 

1.3450 

80.0 

1.0757 

1.0562 

1.0303 

85.0 

.7970 

.7837 

.7667 

90.0 

.5182 

.5112 

.5031 

95.0 

.2394 

.2387 

.2395 

100.0 

-.0394 

-.0339 

-.0242 


5.6817 

3.2193 

1.1727 

.7757 

.3786 

.4481 

.8303 

1.2126 

1.6121 

1.8049 

1.9976 

1.9613 

1.8003 

1.6392 

1.3391 

1.0505 

.7620 

.4921 

.2400 

-.0122 


Table 4.14 




5- 


2 . \ 

n g 3 (n) 

0.0 .0000 

5.0 -.2348 

10.0 -.8710 

15.0 -1.8063 

20.0 -2.9383 

25.0 -4.1649 

30.0 -5.3837 

35.0 -6.4982 

40.0 -7.5155 

45.0 -8.4647 

50.0 -9.3609 

55.0 -10.2188 

60.0 -11.0535 

65.0 -11.8797 

70.0 -12.7146 

75.0 -13.5719 

80.0 -14.4502 

85.0 -15.3469 

90.0 -16.2595 

95.0 -17.1858 

100.0 -18.1231 


gfr) 

85 (n) 

87 (n) 

.0000 

.0000 

.0000 

-.2576 

-.2673 

-.2745 

-.9445 

-.9727 

-.9903 

-1.9320 

-1.9715 

-1.9854 

-3.0915 

-3.1190 

-3.1171 

-4.2941 

-4.2940 

-4.3031 

-5.4478 

-5.4540 

-5.4733 

-6.5450 

-6.5720 

-6.5842 

-7.5840 

-7.6210 

-7.6289 

-8.5629 

-8.5918 

-8.6032 

-9.4797 

-9.5044 

-9.5179 

-10.3493 

-10.3761 

-10.3910 

-11.1963 

-11.2240 

-11.2407 

-12.0341 

-12.0666 

-12.0838 

-12.8759 

-12.9158 

-12.9327 

-13.7350 

-13.7790 

-13.7988 

-14.6201 

-14.6632 

-14.6861 

-15.5278 

-15.5720 

-15.5946 

-16.4533 

-16.5000 

-16.5232 

-17.3917 

-17.4414 

-17.4659 

-18.3385 

-18.3902 

-18.4158 


Table 4.15 
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Figure 4.16 



a 1 



1 00 . 0 


**"*•-_ 

"■S! 

1 

100.0 


Figure 4.17 
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5. Concluding Remarks 

We have presented an approximation theory for numerical solution of the 
discrete-time optimal linear regulator problem in Hilbert space, on both 
finite and infinite time intervals. The motivation for this theory comes from 
optimal control problems for systems involving diffusion equations, hereditary 
differential equations and distributed models of flexible structures. We have 
demonstrated the application of the theory to examples from all three areas 

The solution to the infinite dimensional optimal control problem is based 
on an infinite dimensional Riccati operator equation — a difference equation 
in the finite-time problem and an algebraic equation in the infinite-time 
problem. We have shown that the solution to the infinite dimensional problem 
can be approximated by the solutions to a sequence of finite dimensional 
problems each of which involves a finite dimensional Riccati matrix equation 
to be solved numerically. The finite dimensional problems are just the 
corresponding optimal control problems for finite element approximations to 
the infinite dimensional control system. For the infinite-time problem, the 
finite dimensional Riccati equations usually are solved via eigenspace 
decomposition of the Hamiltonian matrix. 

In both continuous and discrete-time optimal regulator problems for 
distributed systems, the numerical solution often involves solution of large 
Riccati matrix equations. As we observed at the beginning of Section 4, the 
asymptotic relationship between the eigenvalues of a continuous-time 
Hamiltonian system and the eigenvalues of the corresponding discrete— time 
Hamiltonian system is exponential. This means that the approximating finite 
dimensional discrete-time Riccati equations for a given distributed system 
invariably are not as well conditioned as the corresponding continuous-time 
Riccati equations. Nonetheless, as our examples should illustrate, the 
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numerical solution of such problems is well within the reach of current 
computing. To emphasize this, we obtained all of the numerical results in 
this paper on an IBM Personal Computer (not an XT or AT) with 640K of random 
access memory and an Intel 8087 math coprocessor chip. The largest Riccati 
matrix equation that we solved here was a 30 x 30 steady state equation for 
the hub-beam-tip mass example. This solution takes 15 to 20 minutes on the 
PC. We have solved much larger Riccati equations easily on larger mainframe 


computers. 
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